Enhancement of the stabihty of genetic switches by 
overlapping upstream regulatory domains 
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We study genetic switches formed from pairs of mutually repressing operons. The switch stabihty 
is characterised by a well defined lifetime which grows sub-exponentially with the number of copies of 
the most-expressed transcription factor, in the regime accessible by our numerical simulations. The 
stability can be markedly enhanced by a suitable choice of overlap between the upstream regulatory 
domains. Our results suggest that robustness against biochemical noise can provide a selection 
pressure that drives operons, that regulate each other, together in the course of evolution. 
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Biochemical networks are the analog computers of life. 
They allow Uving cells, to detect, transmit, and amplify 
environmental signals, as well as integrate different sig- 
nals in order to recognize patterns in, say, the food sup- 
ply. Indeed, biochemical networks can perform a variety 
of computational tasks analogous to electronic circuits. 
However, their design principles are markedly different. 
In a biochemical network, computations are performed 
by molecules that chemically and physically interact with 
each other. These interactions are stochastic in nature. 
This becomes particularly important when the concen- 
trations are low. In gene regulatory networks, this is 
generally the case: not only the DNA, but also the pro- 
teins that regulate gene expression are often present in 
very small numbers, which can be as low as ten, or even 
fewer. Hence, one would expect that gene regulatory 
networks, in contrast to electronic circuits, are highly 
stochastic and error prone 0, S S • An important 
question, therefore, is how the ability to resist biochemi- 
cal noise constrains the design of the network 0, Q . 

In prokaryotes, the expression of operons — groups of 
contiguous genes that are transcribed into single mRNA 
molecules — is regulated by the binding of transcription 
factors (TPs) to upstream regulatory domains on the 
DNA. A spatial arrangement in which two operons are 
transcribed in diverging directions allows the upstream 
regulatory domains to interfere with each other. This 
affords additional regulatory control. In particular, bio- 
chemical noise in the expression of operons can become 
correlated or anticorrelated. Just as the existence of 
operons provides for correlated gene expression, inter- 
ference between the regulatory domains of two diverging 
operons allows for correlated or anticorrelated expression 
of operons. Here, we show that this can have a dramatic 
influence on the stability of gene regulatory networks. 

Recently, we have performed a statistical analysis of 
the spatial distribution of operons on the genome of E.s- 
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FIG. 1: (a) In a prokaryote, a so-called toggle switch can 
be formed from a pair of operons that mutually repress each 
other, (b) If these operons are transcribed in diverging orien- 
tations, the upstream regulatory domains can interfere with 
each other. 



cherichia coli . The analysis identified a large number 
of motifs in which the operator regions of two operons 
overlap and interfere (6j. Among them are well-known 
examples such as the lysA-lysR and the araBAD-araC 
operon pairs 0- But perhaps the best-known and ar- 
guably the most studied example of such a motif is pro- 
vided by the A-phage switch, which consists of two adja- 
cent operons that mutually repress each other |^. Here, 
we study a minimalist model of such a switch, as shown in 
Fig. ^ In particular, we compare the stability of an 'ex- 
clusive' (XOR) switch, for which the simultaneous bind- 
ing of the repressive TPs for both operons is inhibited, to 
that of a general switch. We find that the exclusive switch 
is much more stable than the general switch. The basis 
reason is that one cannot simultaneously turn off both 
operons in the exclusive switch. This demonstrates the 
potential importance of such motifs in making gene reg- 
ulatory networks robust against biochemical noise. Al- 
though these constructions seem peculiar to prokaryotes, 
evidence for correlated and anticorrelated gene expres- 
sion has also been reported for eukaryotes ^ . 
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TABLE I: Distinct possibilities for the subsets of allowed 
genome states for our switch model. 



case / genome states 
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OA„ 


OB™ 


OA„B 


general 
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exclusive 
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partially co-operative 
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X 
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totally co-operative 


/ 


X 


X 


/ 
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0.6- 



CQ 



ca 
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0.2- 



n = m = 2 




The starting point of our analysis is a set of chemical 
reactions that constitute the switches shown in Fig.^ As 
chemical species, we introduce a pair of TFs which can 
exist as monomers, A and B, or multimers, A„ and Bm. 
The state of the genome is represented by O, 0A„, etc, 
depending on the binding of the TP multimers. Adopting 
a condensed notation in which '|' indicates alternative 
sets of reactants and indicates that the reactants 
are not destroyed by the reaction, the set of chemical 
reactions are 

nA ^ A„, mB ^ B^, {kf, kb) (la) 

O + A„ ^ 0A„, O + B,„ ^ 0B„„ (/con, fcoff) (lb) 
0A„ + B„ I 0B„ + A„ ^ OA„B™, (kon, Ks) (Ic) 
0|OA„-^A, 0|OB™-^B, (A;A),(fcB) (Id) 
A|B^0. (ma),(mb) (le) 

These reactions account for, respectively, the formation 
of multimers, the binding of TF multimers to the genome 
(Eqs. Ijlb|l and fTc| ). the expression of TF monomers, and 
the degradation of TF monomers. Repression of gene ex- 
pression is implicit in Eqs. (jldjl . thus A is expressed if 
and only if B^ is not bound, etc. Reaction rates are as 
indicated, and we define equilibrium constants for mul- 
timerisation, = kf/k^, and binding to the genome, 

Kh = fcon/fcoff- 

Whilst detailed and biologically faithful models can 
be constructed as has been done for the A-phage switch 
[ill H^. the above model is intentionally 'as simple as 
possible'. We believe such an approach is as important 
as detailed biological modelling in elucidating the basic 
physics behind switches. Thus, for example, we have con- 
densed the details of transcription and translation into a 
single reaction step in Eqs. Ijld|l . governed by rate coef- 
ficients and fee- On the other hand, as Cherry and 
Adler have shown the TF binding isotherms must 
satisfy certain criteria in order to make a working switch. 
In the present model this is effected by introducing co- 
operativity through the binding of TF multimers rather 
than monomers. 

In our model the genome is in one of four states 
{0, 0A„, OBm, OA„Bm}- We now include the effect of 
interference between the upstream regulatory domains by 
disallowing some of these states (this is in the spirit of 
simplicity, strictly speaking the effect is to modify the 




FIG. 2: In mean field theory, switching behaviour is confined 
to a wedge in the (/^A/fcA, AtB/^B) plane. Shown are the phase 
diagrams for the dimerising general (solid line) and exclusive 
(dashed line) switches. It is seen that the region of bi-stability 
is larger for the exclusive switch than for the general switch. 
For the partially co-operative dimerising switch (OB2 disal- 
lowed), the wedge moves to //a/^a < 0.10 and fiB/ku < 0.019. 



probabilities of the states). Since the empty genome is 
always a possibility and both A„ and B^ should be al- 
lowed to bind otherwise they would not be TFs, it turns 
out that there are only five possibilities, two of which are 
related by symmetry. The four distinct cases are shown 
in Table and are implemented by excluding some of 
the reactions in Eqs. (|lb|) and l(Tc|) . For example, the 
exclusive switch is obtained by discarding the reactions 
in Eqs. (fTc|l thereby removing the state OA„Bm- 

We first use mean-field theory to analyse the behaviour 
of Eqs. Q. Switching behaviour corresponds to the ap- 
pearance of two distinct stable states in the space of 
TF molecule numbers. Previously, general switches were 
studied in detail by Cherry and Adler , and a specific 
example of the exclusive switch was studied by Kepler 
and Elston fl4|. We extend the analysis of Cherry and 
Adler to determine where switching behaviour can occur, 
for all the cases in Table n Firstly, for n = m = 1, no 
switching behaviour can be found for any case. This con- 
firms that some form of co-operative binding is required. 
For the totally co-operative switch though, switching be- 
haviour cannot be found for any values of n and m. For 
the remaining cases, we have analysed in detail the sit- 
uation for n = m = 2 where both TFs bind as dimers. 
Fig. El shows the regions in the (/^A/fcA, MB/fce) plane 
where switching behaviour is found. Clearly, switching 
behaviour is more extensive for the exclusive switch com- 
pared to the general switch, and is strongly suppressed 
for the partially co-operative switch. Thus we conclude 
that, at least in mean field theory, the structure of the 
switch has a strong influence on the extent of switching 
behaviour. 

To go beyond mean field theory, we have simulated 
the reactions in Eqs. using Gillespie's kinetic Monte- 
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FIG. 3: (a) Typical numbers of TFs as a function of time, (b) 
Cumulative distribution functions for the time intervals be- 
tween zero crossings of A'^a — A^'b . Results are for the dimeris- 
ing exclusive switch at fi/k = 0.45 unless stated otherwise. 



Carlo scheme which generates sample traiectories ap- 
propriate to the chemical master equation [15|. We fo- 
cus on dimerising (n = m = 2) general and exclusive 
switches, and on the symmetry line fcA = fce = fc and 
fiA — MB = We will use the expression rate k « 0.1- 
1 s^^ 1131 s-s a unit of (inverse) time and the degradation 
rate fj, as the main control parameter. The choice of the 
rate constants is biologically motivated, in particular we 
expect expression to be a slow step and the binding equi- 
librium to be biased in favour of bound states [l^l • For 
a baseline set we use k{/V — kiy — ko^/V — Sfcog — 5fc 
{Kd/V = 1 and Kx,/V = 5), where V « 2 ^it? is the 
cell volume [l^ : we assume one copy of the genome is 
present. 

We monitor the total numbers of the TFs, Np^ and A^b, 
including those in dimers and those bound to the genome. 
If the system is behaving as a switch then we typically 
see that one of the TFs is strongly repressed compared to 
the other one. A switching event occurs when the roles 
of the two TFs flip spontaneously, as shown in Fig. 13 

We can obtain more insight into the switching 
behaviour by sampling the probability distribution 
P{Na,Nb) for states in the {Na,Nb) plane, as shown 
in Fig. 21 Switching behaviour appears as a double max- 
imum in probability in this representation, and the tran- 
sition state is seen to lie at low numbers of both TFs. 
Three points are worthy of note. First of all, it is seen 
that the positions of the two stable steady states do not 
depend much on the architecture of the switch. This is 
not surprising, because if one species dominates, both 
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FIG. 4: Probability density in {Na,Nb) plane constructed 
from 2.5 x 10® samples (total duration oi kt ^ 5x 10®), for (a) 
general and (b) exclusive dimerising switches, at fi/k = 0.45. 
Greyscale indicates bin count, logarithmically, from < 1 
(white) to > 10^ (black). Insets show probability density col- 
lapsed onto the A^a — A^b line, plotted as a dimensionless "free 
energy" — log[P{NA — Nb)] (the ordinate zero is arbitrary). 



switches will behave similarly. What is perhaps more 
surprising, is that the pathways for switching are differ- 
ent. The transition paths of the exclusive switch cross 
the transition state surface at higher values of A^a = .^b , 
as compared to the general switch. The reason is that in 
the general switch both genes can be repressed simulta- 
neously, while in the exclusive switch only one gene can 
be turned off at the time. More importantly, however, 
the barrier for flipping the switch is higher for the exclu- 
sive switch than for the general switch, as can be seen in 
the insets in Fig. ^ This is because for a switch to flip, 
two events have to happen. First of all, the system has to 
wait for a rare fluctuation by which the concentration of 
the dominant species decreases; this allows for the syn- 
thesis of the other component. Subsequently, the latter 
component has to bind to its operator site in order to 
toggle the switch. In the general switch, the latter event 
is more probable, because the minor component can bind 
to its site as soon as it is synthesized, while in the ex- 
clusive switch the dominant species first has to dissociate 
from the DNA. This is the main reason why the exclusive 
switch is more stable than the general switch. 

We have also characterised the switching dynamics, by 
constructing the cumulative distribution function F{At) 
for the time intervals At between zero crossings of the 
order parameter Na — Nb ■ About 50% of F arises from 
noise on a time scale At ~ k"^ as the system jitters 
around the transition state, but for At ^ k~^, and pro- 
vided we are well into the switching regime, we invari- 
ably see Poisson statistics with F ^ 1 — exp[— Ai/r] (see 
Fig. Etb)). This firstly confirms that the switch states 
have a well defined lifetime r, and secondly allows us to 
extract an accurate estimate of the value of r. 

Bialek has suggested that the switch lifetime may grow 
exponentially with the number of molecules involved in 
switching between states Motivated by this, we 

monitor the mean number N of the most-expressed TF, 
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FIG. 5: Switch lifetime as a function of mean number of 
most-expressed TF, for the general (solid line) and exclusive 
(dashed line) cases. The exclusive switch becomes orders of 
magnitude more stable than the general switch at high num- 
bers of the expressed TF. 



simulations suggest that robustness against biochemical 
noise can provide a selection pressure that drives pairs of 
operons, that either regulate each other or are controlled 
by a common transcription factor, towards each other in 
the course of evolution. 
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defined to be the time average of max(A^A, ^b)- We can 
also calculate N from mean field theory, and we find good 
agreement between this and the value measured in the 
simulations, as /i varies. 

Qualitative support for Bialek's conjecture comes from 
Fig.El It shows that r grows very rapidly with N, which 
is the basic reason why extremely stable switches can be 
built with at most a few hundred expressed proteins. In 
contrast to Bialek's conjecture, however, t{N) appears to 
be sub-exponential in N. Interestingly, we can fit t{N) to 
a form that is suggested by an analysis of a related prob- 
lem, namely that of switching between broken-symmetry 
phases in a driven diffusive model 0, 0| . This suggests 
that the ultimate scaling is r ~ iV" exp[67V], where a and 
b are constants. Note that this corresponds to Bialek's 
conjecture, but with a logarithmic correction in N. 

More importantly, however. Fig. |5l clearly demon- 
strates that the switch construction does indeed have a 
marked influence on the stability of the switch. It shows 
that the lifetime of the exclusive switch grows much more 
rapidly with mean copy number than that of the general 
switch Our simulations cover 10 < TV < 30, but if 
we extrapolate our results to TV « 100, then kr w 10^- 
10^ for the general switch but fcr « 10^-10^" for the 
exclusive switch. In the latter case, this corresponds to 
lifetimes measured in tens of years. Such extremely long 
lifetimes have been reported for A-phage 12]. 

In summary, a genetic switch is intrinsically stochas- 
tic, because of the molecular character of its components. 
However, our simulations demonstrate that the stability 
of a genetic switch can be strongly enhanced by spatially 
arranging the operons such that competing regulatory 
molecules mutually exclude each other at the operator 
regions. Such a spatial arrangement can be achieved if 
the two operons lie next to each other on the DNA and 
are transcribed in diverging directions - a network mo- 
tif that has been identified by our statistical analysis of 
the gene regulatory network of E. coli 6]. Hence, our 
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